#include<iostream>
#include<fstream>
#include<cmath>
#include<string>
#include<chrono>
#include <omp.h>
using namespace std;
using namespace chrono;
#define N 8399 //矩阵列数
#define C2 4535 //被消元行行数
int s[N][N] = { 0 }, a[C2][N] = { 0 }, reach[C2] = {0};
int p = 5; //每次导入的消元子数
int thread_num = 8, chunk_size = 1;

void initialize()
{
	for (int i = 0; i < N; i++)
	{
		if (i < C2)
			reach[i] = 0;
		for (int j = 0; j < N; j++)
		{
			if (i < C2)
				a[i][j] = 0;
			s[i][j] = 0;
		}
	}
	int t1, t2,i;
	fstream fin1("xyz.txt");
	fin1 >> t1;
	i = N - t1 - 1;
	s[i][i] |= 1;
	while (!fin1.eof())
	{
		fin1 >> t2;
		if (t2 > t1)
		{
			i = N - t2 - 1;
			//s[i][i] = 1;
			s[i][i] |= 1;
			t1 = t2;
		}
		else
		{
			s[i][N - t2 - 1] |= 1;
			t1 = t2;
		}
	}
	fin1.close();

	fstream fin2("bxyh.txt");
	i = 0;
	fin2 >> t1;
	a[i][N - t1 - 1] |= 1;
	while (!fin2.eof())
	{
		fin2 >> t2;
		if (t2 > t1)
		{
			a[++i][N - t2 - 1] |= 1;
			t1 = t2;
		}
		else
		{
			a[i][N - t2 - 1] |= 1;
			t1 = t2;
		}
	}
	fin2.close();
}



void print()
{
	for (int i = 0; i <C2; i++)
	{
		for (int j = 0; j < N; j++)
			cout << a[i][j] << " ";
		cout << endl;
	}
	cout << endl;
}

void LU()
{
	int i, j, k, m;
	for (i = N - 1; i >= 0; i -= p)
	{
		if (i - p < 0)
			p = i + 1;
		for (j = 0; j < C2; j++)
		{
			if (reach[j] == 0)
				for (k = 0; k < i + p; k++)
				{
					if (a[j][k] == 1 && s[k][k] == 1)
					{
						for (m = k; m < N; m++)
						{
							a[j][m] = a[j][m] ^ s[k][m];
						}
					}
					else if (a[j][k] == 1 && s[k][k] == 0)
					{
						for (m = k; m < N; m++)
						{
							s[k][m] = a[j][m];
						}
						reach[j] = 1;
						break;
					}
				}
		}
	}
}

void LU_OMP()
{
	int i, j, k, m;
#pragma omp parallel num_threads(thread_num), private(i, j, k,m)
	for (i = N -1; i>=0; i -= p)
	{
		if (i - p < 0) 
			p = i + 1;
#pragma omp for schedule(static,chunksize)
		for (j = 0; j < C2; j++)
		{
			if(reach[j]==0)
			for (k = 0; k < i+p; k++)
			{
				if (a[j][k] == 1 &&s[k][k]==1)
				{
					for (m = k; m < N; m++)
					{
						a[j][m] = a[j][m] ^ s[k][m];
					}
				}
				else if (a[j][k] == 1 && s[k][k] == 0)
				{
					for (m = k; m < N; m++)
					{
						s[k][m] = a[j][m];
					}
					reach[j] = 1;
					break;
				}				
			}
		}
	}
}





void store()
{
	ofstream fout("result.txt");
	for (int i = 0; i < C2;i++)
	{
		bool flag = false;
		for (int j = 0; j < N; j++)
		{
			if (a[i][j] == 1)
			{
				fout << N - j - 1<<" ";
				flag = true;
			}
		}
		if (flag)
			fout << endl;
	}
}


void timing(void(*func)())
{
	initialize();
	typedef std::chrono::high_resolution_clock Clock;
	auto t1 = Clock::now();
	func();
	auto t2 = Clock::now();
	store();
	cout << std::chrono::duration_cast<std::chrono::nanoseconds>(t2 - t1).count() / 1e+6 << endl;
}


int main()
{
	timing(LU);
	timing(LU_OMP);
}